*Sydney computer /Users/kongn/Dropbox (Sydney Uni)/Self Control Stability Project/SOEP/DO
global SOEP "/Users/wkon4921/Dropbox (Sydney Uni)/Self Control Stability Project/SOEP/"
*QUT computer
global SOEP "/Users/kongn/Dropbox (Sydney Uni)/Self Control Stability Project/SOEP/"

use "${SOEP}OUTDATA/inno2017_2020_ppfad_V3",clear /*Long form from V3*/
cd "${SOEP}OUTPUT/"

di "****************************************************"
di "**************Figure 1 Dif_SC Age  *****************"
di "****************************************************"
use "${SOEP}OUTDATA/inno2017_2020_ppfad_V3",clear /*Long form from V3*/
cd "${SOEP}OUTPUT/"

twoway (histogram SC if year==2017,width(2) color(green)) ///
       (histogram SC if year==2020,width(2) fcolor(none) lcolor(black)) ///
	   (kdensity SC if year==2017,lpattern(dash) lwidth(thick) lc(grey)) (kdensity SC if year==2020, lpattern(solid) lwidth(thick) lc(black)) ///
	   ,xtitle("SC index")legend(order(1 "2017" 2 "2020" 3 "Kernel 2017" 4 "Kernel 2020"  )) xsc(r(18 65)) xlabel(18[2]65) 
graph save "${SOEP}OUTPUT/SC2017_2020_V3",replace
graph export "${SOEP}OUTPUT/SC2017_2020_V3.pdf",replace  
	   twoway (histogram SC if year==2017,width(4) color(green)) ///
       (histogram SC if year==2020,width(4) fcolor(none) lcolor(black)) ///
	   (kdensity SC if year==2017,lpattern(dash) lwidth(thick) lc(grey)) (kdensity SC if year==2020, lpattern(solid) lwidth(thick) lc(black)) ///
	   ,xtitle("SC index")legend(order(1 "2017" 2 "2020" 3 "Kernel 2017" 4 "Kernel 2020"  )) xsc(r(18 65)) xlabel(18[4]65) 
graph save "${SOEP}OUTPUT/SC2017_2020_V4",replace
graph export "${SOEP}OUTPUT/SC2017_2020_V4.pdf",replace
	   *xlabel(minmax)
di "****************************************************"
di " FIGURE 2 noise to signal ratio"
di "****************************************************"

use  "${SOEP}OUTDATA/noisetosignal.dta",clear
save "${SOEP}OUTDATA/noisetosignal.dta",replace

/*
gen R=0
rho_0=0.67*(1+R)
*/

twoway (line noisetosignalratio truecorrelation)

gen Rho_0=0.67*(1+within_sigma/among_sigma)
replace  Rho_0=0.67*(1+within_sigma/among_sigma)
replace Rho_0=1 if Rho_0>1

twoway (line TrueRho NSR), xtitle(Noise-to-signal ratio) ytitle(True correlation) yline(0.67)
graph save "${SOEP}/OUTPUT/measurement_error_V2",replace
graph export "${SOEP}/OUTPUT/measurement_error_V2.pdf",replace

di "****************************************************"
di " FIGURE 3 use google chart"
di "****************************************************"
 di "use google chart, webiste https://jsfiddle.net/0x67vrp9/"
 google.charts.load('current', {'packages':['sankey']});
      google.charts.setOnLoadCallback(drawChart);

      function drawChart() {
        var data = new google.visualization.DataTable();
        data.addColumn('string', 'From');
        data.addColumn('string', 'To');
        data.addColumn('number', 'Weight');
        data.addRows([
          [ 'Q1', 'q1', 13.5 ],
          [ 'Q1', 'q2', 5.25 ],
          [ 'Q1', 'q3', 2.83 ],
          [ 'Q1', 'q4', 0.4 ],
          [ 'Q2', 'q1', 5.98 ],
          [ 'Q2', 'q2', 7.52 ],
          [ 'Q2', 'q3', 6.55 ],
          [ 'Q2', 'q4', 2.51 ],        
          [ 'Q3', 'q1', 2.83],
          [ 'Q3', 'q2', 7.03],
          [ 'Q3', 'q3', 11.16],
          [ 'Q3', 'q4', 7.84],
          [ 'Q4', 'q1', 0.97],
          [ 'Q4', 'q2', 2.18],
          [ 'Q4', 'q3', 6.87],
          [ 'Q4', 'q4', 16.59]
        ]);

        // Sets chart options.
        var options = {
          width: 300,
        };

        // Instantiates and draws our chart, passing in some options.
        var chart = new google.visualization.Sankey(document.getElementById('sankey_basic'));
        var colors = ['#cab2d6', '#cab2d6', 
        							'#fb9a99', '#ffff99',
                 			'#1f78b4', '#fb9a99', 
                      '#ffff99', '#1f78b4'];

    var options = {
      height: 400,
      sankey: {
        node: {
          colors: colors
        },
        link: {
          colorMode: 'gradient',
          colors: colors
        }
      }
    };
        chart.draw(data, options);
      }

di "****************************************************"
di "**************Figure 4 Dif_SC Age  *****************"
di "****************************************************"

use "${SOEP}OUTDATA/inno2017_2020_ppfad_V3",clear /*Long form from V3*/

di "use long form directly, adding linear prediction"
	foreach v in dif_SC {
						preserve
						keep if year==2017
						de,s
						collapse (mean) mean`v'= `v' (sd) sd`v'=`v' (count) n=`v', by(agegroup)
						generate hi`v' = mean`v' + invttail(n-1,0.025)*(sd`v' / sqrt(n)) 
						generate lo`v' = mean`v' - invttail(n-1,0.025)*(sd`v' / sqrt(n))
						#delimit ;
						twoway 
						 (rarea hi`v' lo`v' agegroup, fc(lavender%10) sort astyle(ci)) 
						 (line mean`v' agegroup,sort lpattern(dash) lc(navy)) 
						(function y= 0, range(0 11) lc(cyan) lpattern(solid) lwidth(thick)) 
						, xlabel(			
												0 "17-24"  			
												1 "25-29"           
												2 "30-34"           
												3 "35-39"           
												4 "40-44"           
												5 "45-49"           	
												6 "50-54"					
												7 "55-59"              
												8 "60-64"			
												9 "65-69"           
												10 "70-74"					
												11 ">=75"       
																	   ) 	   
						legend(label(1 "95% CI") label(2 "Mean level change in SC") label(3 "y=0")) 
						yline(2.925334 -3.023639) ytitle("Changes in SC index")  ysc(r(-10 10)) 
						ylabel(-10 "-10"  
								-8 "-8"   
								-6 "-6"   
								-4 "-4"   
								-2 "-2"   
								0  "0"    
								2 "2"     
								4 "4"     
								6 "6"     
								8 "8"     
								10 "10");
						#delimit cr
						graph save "`v'_agegroup_y0_V6.gph",replace
						graph export "`v'_agegroup_y0_V6.pdf",replace
						restore
						}

				
di "******************OTHER STATS*******"
di "***********************************"	   				
 log using "${SOEP}OUTPUT/percent change by sd_v2.smcl",replace			

di "% of standard deviation change "" 
di "*****FOR SAMPLE 25-65*******************"	
use "${SOEP}OUTDATA/inno2017_2020_ppfad_V3",clear
keep if year==2017 /*so one obs represents one person */
di "********************"
di "For full age sample"
di "********"
di "% not changed"
tab dif_SC 
di "********"
di "changed less than 1 sd"
sum SC2017 
local sd=r(sd)
local N=r(N)
*di r(sd)
di "-`sd',`sd'"
sum dif_SC if inrange(dif_SC,-`sd',`sd') & dif_SC!=0
local i=r(N)
di `i' " respondents changed within 1 sd"
di `i'/`N'  "%respondents changed within 1 sd"

di "********"
di "changed more than 1 sd less than 2sd"
sum SC2017 
local sd=r(sd)
local N=r(N)
*di r(sd)
sum dif_SC if dif_SC>`sd' & dif_SC<=2*`sd'
local i=r(N)
sum dif_SC if dif_SC<-`sd' & dif_SC>=-2*`sd'
local j=r(N)
di `i'+`j' " respondents changed more than 1 sd"
di (`i'+`j')/`N'  "%respondents changed more than 1 sd"



di "********"
di "changed more than 2 sd less than 3sd"
sum SC2017 
local sd=r(sd)
local N=r(N)
*di r(sd)
sum dif_SC if dif_SC>2*`sd' & dif_SC<=3*`sd'
local i=r(N)
sum dif_SC if dif_SC<-2*`sd' & dif_SC>=-3*`sd'
local j=r(N)
di `i'+`j' " respondents changed more than 2 sd less than 3sd"
di (`i'+`j')/`N'  "%respondents more than 2 sd less than 3sd"

di "********"
di "changed more than 3 sd"
sum SC2017
local sd3=r(sd)*3
local N=r(N)
*di r(sd)*3

sum dif_SC if dif_SC>`sd3'
local i=r(N)
sum dif_SC if dif_SC<-`sd3'
local j=r(N)
di `i'+`j' "respondents changed more than 3 sd"
di (`i'+`j')/`N' " %respondents changed more than 3 sd"




        
log close       
        
        
        
        
        




